#### 03_prepareRollCallBrazil-1.R ####
rm(list=ls())

print (version)

# platform       aarch64-apple-darwin20      
# arch           aarch64                     
# os             darwin20                    
# system         aarch64, darwin20           
# status                                     
# major          4                           
# minor          1.2                         
# year           2021                        
# month          11                          
# day            01                          
# svn rev        81115                       
# language       R                           
# version.string R version 4.1.2 (2021-11-01)
# nickname       Bird Hippie            

# Load the following libraries
library (openxlsx)
library (reshape2)
library (gtools)
library (car)
library (plyr)
library (MCMCpack)
library (zoo)
library (wnominate)
library (pscl)
library (runjags)
library (random)

# This file prepares a long-format dataset
# with information on roll-call votes, 
# vote characteristics, and legislator characteristics.

# Do not change the following seed
set.seed (10882)

# Set working directory
savePathRep <- "/all_data/"
setwd (savePathRep)

base <- read.table ("brazil_original_rollcall.txt", sep=",", header=TRUE, fileEncoding="UTF-8")
year <- as.numeric (substr (base$Id_Votacao, 1, 4))

# Eliminate votes prior to 2007
smallBase <- base[year>=2007,]
rm (base)

# Create dataset for legislature 54
Leg54 <- smallBase[is.element (substr (smallBase$Id_Votacao, 1, 4), 2011:2014),]
rm (smallBase)

# Separate vote-specific information from the rest of the set
vote.characteristics <- c("Id_Votacao","Data","Procedim","Assunto","PROJ_ANO","Projeto","Unanime")

Votacoes.54 <- Leg54[,is.element(colnames(Leg54),vote.characteristics)] # Each vote repeated as many times as there are deputies
Votacoes.54 <- Votacoes.54[duplicated(Votacoes.54)==FALSE,]

# Functions employed below to find party switchers
is.unique <- function (x) {
  l <- length (unique (na.omit(x)))
  un <- ifelse (l==1, 1, 0 )
  return (un)
}

modal.party <- function (x) {
  p <- unique (na.omit(x))[1]
  return (p)
}


########################################################################
# Legislature 54
# Separate legislator-specific information from the rest of the set
########################################################################
leg.characteristics <- c("Id_Votacao"
                         , colnames(Leg54)[!is.element(colnames(Leg54),vote.characteristics)])
Deputados.54 <- Leg54[, is.element(colnames(Leg54), leg.characteristics)]
ind <- sapply(Deputados.54, is.factor) 
Deputados.54[ind] <- lapply(Deputados.54[ind], "[", drop=TRUE)

w.vote.54  <- reshape(Deputados.54[,c(1,3,5)], 
                      timevar = "Id_Votacao",
                      v.names = "Voto",
                      idvar = "Nome_Deputado",
                      direction = "wide")
rownames (w.vote.54) <- 1:nrow(w.vote.54)

w.party.54 <- reshape(Deputados.54[,c(1,3,4)], 
                      timevar = "Id_Votacao",
                      v.names = "Sigla_Partido",
                      idvar = "Nome_Deputado",
                      direction = "wide")
rownames (w.party.54) <- 1:nrow(w.party.54)


## Add deputy information 
deputy.Info <- read.table ("brazil_deputy_information_leg54.txt", header=T, sep="\t")
deputy.Info <- deputy.Info[deputy.Info$RollcallName != "NOT IN LEG",]

is.element(unique(deputy.Info$RollcallName), unique(w.party.54$Nome_Deputado))
is.element(unique(w.party.54$Nome_Deputado), unique(deputy.Info$RollcallName))

# all legislators are in Brazilian deputies dataset
w.party.54$Nome_Deputado[!is.element(w.party.54$Nome_Deputado, deputy.Info$RollcallName)]

# 29 legislators in Brazilian deputies dataset are not in roll-call dataset 
deputy.Info$RollcallName[!is.element(deputy.Info$RollcallName, w.party.54$Nome_Deputado)]


## Convert factors to character strings
w.party.54[] <- lapply(w.party.54, as.character)

same.party.54 <- apply (w.party.54[,-1], 1, is.unique)
avg.party.54  <- apply (w.party.54[,-1], 1, modal.party)
avg.party.54  <- ifelse (avg.party.54=="PFL", "DEM", avg.party.54)
party.colors  <- car::recode (avg.party.54, " 'DEM'='blue'
                         ; 'PMDB'='purple'; 'PSDB'='black'
                         ; 'PT'='red'; 'PP'='orange'
                         ; 'PR'='yellow'; else='gray'")

# repeat.NA will fill in missing party labels with the last observed party label
# repeat.NA function from material downloaded from
# https://stackoverflow.com/questions/7735647/replacing-nas-with-latest-non-na-value#7735681
# on June 22, 2018
replace.NA <- function (tmp) {
  z  <- !is.na(tmp)       # indicates the positions of y whose values we do not want to overwrite
  z  <- z | !cumsum(z)    # for leading NA's in y, z will be TRUE, otherwise it will be FALSE where y has a NA and TRUE where y does not have a NA
  y  <- tmp[z][cumsum(z)]
  return (y)
}

# filled up party matrix
party.only <- w.party.54[,-1]
full.w.party.54 <- t (apply (party.only, 1, replace.NA))

# Build from the following to eliminate "party labels" if 
# deputy has left congress
# we eliminate deputies that are absent more than 10 votes at the end of 
# the congressional term

NA.if.legislator.out <- function (x, y) {
  NonNAindex <- which(!is.na(x))
  lastNonNA <- max(NonNAindex)
  if (lastNonNA < (ncol(party.only)-10))  {
    y[(lastNonNA+1):length(y)] <- NA
  } 
  return (y)
}

correct.w.party.54 <- matrix (NA, ncol=ncol(full.w.party.54), nrow=nrow(full.w.party.54))
for (i in 1:nrow(full.w.party.54)) {
  correct.w.party.54[i,] <- NA.if.legislator.out(x=party.only[i,], y=full.w.party.54[i,])
}
colnames (correct.w.party.54) <- colnames (w.party.54)[-1]

#### Build RollCall matrix, substituting numbers for letters ####
f <- apply (w.vote.54, 2, plyr::revalue, c("A"=-9, "F"=-7, "O"=-5, "N"=0, "S"=1))
simple.rc <- apply (f, 2, as.numeric)
rm (f)
rc <- simple.rc[,-grep("Nome_Deputado", colnames(simple.rc))]


f2 <- apply (w.vote.54, 2, plyr::revalue, c("A"=NA, "F"=NA, "O"=2, "N"=0, "S"=1))
rc.s <- apply (f2, 2, as.numeric)
rm (f2)
rc.2 <- rc.s[,-grep("Nome_Deputado", colnames(rc.s))]

# Extract vote information
propAye  <- colMeans (rc.s[,-1], na.rm=T)
nPresent <- colSums (!is.na(rc.s[,-1]))
vote.characteristics <- c("Id_Votacao","Data","Procedim","Assunto","PROJ_ANO","Projeto","Unanime")
Votes.54 <- Leg54[,is.element(colnames(Leg54), vote.characteristics)]
Votes.54 <- Votes.54[!duplicated(Votes.54),]
Votes.54$Data <- gsub (' 0:00:00', '', x=Votes.54$Data)
Votes.54$Unanime <- ifelse (Votes.54$Unanime=="S", 1, 0)
Votes.54$propAye <- round (propAye, 2)
Votes.54$present <- nPresent

# Read additional information for votes in Leg 54
session.info <- read.table (file="brazil_vote_information_leg54.txt"
                            , sep="\t", header=T)
identical (Votes.54$Id_Votacao, as.numeric (session.info$Id_Votacao))
Votes.54$session <- session.info$Session
rm (session.info)

parties <- c("PFL","PT","PSDB","PMDB","PTB","PDT","PP","PPS","PR","PSB","PV1","PCdoB","DEM")


# Import party pressure data for PT, to find out IDs of bills we need 

vote.code<- read.table("brazil_votes_codification.txt")
valid.vote <- vote.code$vote.number
which.vote <- paste0 ("Voto.", valid.vote)
which.party <- paste0 ("Sigla_Partido.", valid.vote)
#rm (dataPT)
rm (vote.code)

unityRollCall <- rc.2[,is.element (colnames (rc.2), which.vote)]
unityParty <- correct.w.party.54[,is.element (colnames (correct.w.party.54), which.party)]
lopsided.votes <- abs(propAye[is.element(colnames (rc.2), which.vote)] - 0.5) + 0.5


aye.vote.dem <- aye.vote.pt <- aye.vote.pr <- aye.vote.pp <- aye.vote.psdb <- aye.vote.pmdb <- c()
nay.vote.dem <- nay.vote.pt <- nay.vote.pr <- nay.vote.pp <- nay.vote.psdb <- nay.vote.pmdb <- c()
abs.vote.dem <- abs.vote.pt <- abs.vote.pr <- abs.vote.pp <- abs.vote.psdb <- abs.vote.pmdb <- c()
obs.vote.dem <- obs.vote.pt <- obs.vote.pr <- obs.vote.pp <- obs.vote.psdb <- obs.vote.pmdb <- c()

for (i in 1:ncol(unityParty)){
  aye.vote.dem[i] <- length(unityParty[,i][unityParty[,i]=="DEM" & !is.na(unityParty[,i]) & unityRollCall[,i]==1 & !is.na(unityRollCall[,i])])
  aye.vote.pt[i] <- length(unityParty[,i][unityParty[,i]=="PT" & !is.na(unityParty[,i]) & unityRollCall[,i]==1 & !is.na(unityRollCall[,i])])
  aye.vote.pr[i] <- length(unityParty[,i][unityParty[,i]=="PR" & !is.na(unityParty[,i]) & unityRollCall[,i]==1 & !is.na(unityRollCall[,i])])
  aye.vote.pp[i] <- length(unityParty[,i][unityParty[,i]=="PP" & !is.na(unityParty[,i]) & unityRollCall[,i]==1 & !is.na(unityRollCall[,i])])
  aye.vote.psdb[i] <- length(unityParty[,i][unityParty[,i]=="PSDB" & !is.na(unityParty[,i]) & unityRollCall[,i]==1 & !is.na(unityRollCall[,i])])
  aye.vote.pmdb[i] <- length(unityParty[,i][unityParty[,i]=="PMDB" & !is.na(unityParty[,i]) & unityRollCall[,i]==1 & !is.na(unityRollCall[,i])])
}

for (i in 1:ncol(unityParty)){
  nay.vote.dem[i] <- length(unityParty[,i][unityParty[,i]=="DEM" & !is.na(unityParty[,i]) & unityRollCall[,i]==0 & !is.na(unityRollCall[,i])])
  nay.vote.pt[i] <- length(unityParty[,i][unityParty[,i]=="PT" & !is.na(unityParty[,i]) & unityRollCall[,i]==0 & !is.na(unityRollCall[,i])])
  nay.vote.pr[i] <- length(unityParty[,i][unityParty[,i]=="PR" & !is.na(unityParty[,i]) & unityRollCall[,i]==0 & !is.na(unityRollCall[,i])])
  nay.vote.pp[i] <- length(unityParty[,i][unityParty[,i]=="PP" & !is.na(unityParty[,i]) & unityRollCall[,i]==0 & !is.na(unityRollCall[,i])])
  nay.vote.psdb[i] <- length(unityParty[,i][unityParty[,i]=="PSDB" & !is.na(unityParty[,i]) & unityRollCall[,i]==0 & !is.na(unityRollCall[,i])])
  nay.vote.pmdb[i] <- length(unityParty[,i][unityParty[,i]=="PMDB" & !is.na(unityParty[,i]) & unityRollCall[,i]==0 & !is.na(unityRollCall[,i])])
}

for (i in 1:ncol(unityParty)){
  obs.vote.dem[i] <- length(unityParty[,i][unityParty[,i]=="DEM" & !is.na(unityParty[,i]) & unityRollCall[,i]==2 & !is.na(unityRollCall[,i])])
  obs.vote.pt[i] <- length(unityParty[,i][unityParty[,i]=="PT" & !is.na(unityParty[,i]) & unityRollCall[,i]==2 & !is.na(unityRollCall[,i])])
  obs.vote.pr[i] <- length(unityParty[,i][unityParty[,i]=="PR" & !is.na(unityParty[,i]) & unityRollCall[,i]==2 & !is.na(unityRollCall[,i])])
  obs.vote.pp[i] <- length(unityParty[,i][unityParty[,i]=="PP" & !is.na(unityParty[,i]) & unityRollCall[,i]==2 & !is.na(unityRollCall[,i])])
  obs.vote.psdb[i] <- length(unityParty[,i][unityParty[,i]=="PSDB" & !is.na(unityParty[,i]) & unityRollCall[,i]==2 & !is.na(unityRollCall[,i])])
  obs.vote.pmdb[i] <- length(unityParty[,i][unityParty[,i]=="PMDB" & !is.na(unityParty[,i]) & unityRollCall[,i]==2 & !is.na(unityRollCall[,i])])
}

for (i in 1:ncol(unityParty)){
  abs.vote.dem[i] <- length(unityParty[,i][unityParty[,i]=="DEM" & !is.na(unityParty[,i]) & is.na(unityRollCall[,i])])
  abs.vote.pt[i] <- length(unityParty[,i][unityParty[,i]=="PT" & !is.na(unityParty[,i]) & is.na(unityRollCall[,i])])
  abs.vote.pr[i] <- length(unityParty[,i][unityParty[,i]=="PR" & !is.na(unityParty[,i]) & is.na(unityRollCall[,i])])
  abs.vote.pp[i] <- length(unityParty[,i][unityParty[,i]=="PP" & !is.na(unityParty[,i]) & is.na(unityRollCall[,i])])
  abs.vote.psdb[i] <- length(unityParty[,i][unityParty[,i]=="PSDB" & !is.na(unityParty[,i]) & is.na(unityRollCall[,i])])
  abs.vote.pmdb[i] <- length(unityParty[,i][unityParty[,i]=="PMDB" & !is.na(unityParty[,i]) & is.na(unityRollCall[,i])])
}


prop.abs.dem <- prop.abs.pt <- prop.abs.pr <- prop.abs.pp <- prop.abs.psdb <- prop.abs.pmdb <- c()
prop.aye.dem <- prop.aye.pt <- prop.aye.pr <- prop.aye.pp <- prop.aye.psdb <- prop.aye.pmdb <- c()
prop.aye.dem.no.abs <- prop.aye.pt.no.abs <- prop.aye.pr.no.abs <- prop.aye.pp.no.abs <- prop.aye.psdb.no.abs <- prop.aye.pmdb.no.abs <- c()

for (i in 1:ncol(unityParty)){
  prop.abs.dem <- abs.vote.dem / rowSums (cbind (aye.vote.dem, nay.vote.dem, abs.vote.dem, obs.vote.dem))
  prop.abs.pt <- abs.vote.pt / rowSums (cbind (aye.vote.pt, nay.vote.pt, abs.vote.pt, obs.vote.pt))
  prop.abs.pr <- abs.vote.pr / rowSums (cbind (aye.vote.pr, nay.vote.pr, abs.vote.pr, obs.vote.pr))
  prop.abs.pp <- abs.vote.pp / rowSums (cbind (aye.vote.pp, nay.vote.pp, abs.vote.pp, obs.vote.pp))
  prop.abs.psdb <- abs.vote.psdb / rowSums (cbind (aye.vote.psdb, nay.vote.psdb, abs.vote.psdb, obs.vote.psdb))
  prop.abs.pmdb <- abs.vote.pmdb / rowSums (cbind (aye.vote.pmdb, nay.vote.pmdb, abs.vote.pmdb, obs.vote.pmdb))
}

for (i in 1:ncol(unityParty)){
  prop.aye.dem <- aye.vote.dem / rowSums (cbind (aye.vote.dem, nay.vote.dem, abs.vote.dem, obs.vote.dem))
  prop.aye.pt <- aye.vote.pt / rowSums (cbind (aye.vote.pt, nay.vote.pt, abs.vote.pt, obs.vote.pt))
  prop.aye.pr <- aye.vote.pr / rowSums (cbind (aye.vote.pr, nay.vote.pr, abs.vote.pr, obs.vote.pr))
  prop.aye.pp <- aye.vote.pp / rowSums (cbind (aye.vote.pp, nay.vote.pp, abs.vote.pp, obs.vote.pp))
  prop.aye.psdb <- aye.vote.psdb / rowSums (cbind (aye.vote.psdb, nay.vote.psdb, abs.vote.psdb, obs.vote.psdb))
  prop.aye.pmdb <- aye.vote.pmdb / rowSums (cbind (aye.vote.pmdb, nay.vote.pmdb, abs.vote.pmdb, obs.vote.pmdb))
}

for (i in 1:ncol(unityParty)){
  prop.aye.dem.no.abs <- aye.vote.dem / rowSums (cbind (aye.vote.dem, nay.vote.dem))
  prop.aye.pt.no.abs <- aye.vote.pt / rowSums (cbind (aye.vote.pt, nay.vote.pt))
  prop.aye.pr.no.abs <- aye.vote.pr / rowSums (cbind (aye.vote.pr, nay.vote.pr))
  prop.aye.pp.no.abs <- aye.vote.pp / rowSums (cbind (aye.vote.pp, nay.vote.pp))
  prop.aye.psdb.no.abs <- aye.vote.psdb / rowSums (cbind (aye.vote.psdb, nay.vote.psdb))
  prop.aye.pmdb.no.abs <- aye.vote.pmdb / rowSums (cbind (aye.vote.pmdb, nay.vote.pmdb))
}

# Votes with very high rates of abstention among some parties
dem.all.abs <- colnames(rc)[prop.abs.dem>0.7]
psdb.all.abs <- colnames(rc)[prop.abs.psdb>0.7]
pt.all.abs <- colnames(rc)[prop.abs.pt>0.7]
pmdb.all.abs <- colnames(rc)[prop.abs.pmdb>0.7]
pr.all.abs <- colnames(rc)[prop.abs.pr>0.7]
pp.all.abs <- colnames(rc)[prop.abs.pp>0.7]

# Find votes where government parties vote together against opposition parties
all.no.abs <- data.frame (cbind(seq(1,length(prop.aye.dem),by=1), prop.aye.pmdb.no.abs, prop.aye.pt.no.abs, prop.aye.pp.no.abs, prop.aye.dem.no.abs, prop.aye.pr.no.abs, prop.aye.psdb.no.abs))
all.no.abs$Vote.Number <- paste0 ("Voto.", substr(colnames(unityParty), start=15, stop=nchar(colnames(unityParty))))
all.no.abs <- na.omit (all.no.abs)

all.no.abs[all.no.abs$Vote.Number=="Voto.2011043",]

t<-.8      #threshold
# Government in favor, opposition against
round (all.no.abs[all.no.abs$prop.aye.pmdb.no.abs > t & all.no.abs$prop.aye.pt.no.abs > t & all.no.abs$prop.aye.pp.no.abs > t & all.no.abs$prop.aye.dem.no.abs < (1-t) & all.no.abs$prop.aye.psdb.no.abs < (1-t) & all.no.abs$prop.aye.pr.no.abs > t, -ncol(all.no.abs)], 1)

# Government against, opposition in favor
round (all.no.abs[all.no.abs$prop.aye.pmdb.no.abs < (1-t) & all.no.abs$prop.aye.pt.no.abs < (1-t) & all.no.abs$prop.aye.pp.no.abs < (1-t) & all.no.abs$prop.aye.dem.no.abs > t & all.no.abs$prop.aye.psdb.no.abs > t & all.no.abs$prop.aye.pr.no.abs < (1-t), -ncol(all.no.abs)], 1)

# Left-center together (PT, PMDB, PSDB), against right (DEM, PR, PP)
round (all.no.abs[all.no.abs$prop.aye.pmdb.no.abs >t & all.no.abs$prop.aye.pt.no.abs > t & all.no.abs$prop.aye.psdb.no.abs > t & all.no.abs$prop.aye.dem.no.abs < (1-t) & all.no.abs$prop.aye.pp.no.abs < (1-t) & all.no.abs$prop.aye.pr.no.abs < (1-t), -ncol(all.no.abs)], 1)

round (all.no.abs[all.no.abs$prop.aye.pmdb.no.abs < (1-t) & all.no.abs$prop.aye.pt.no.abs < (1-t) & all.no.abs$prop.aye.psdb.no.abs < (1-t) & all.no.abs$prop.aye.dem.no.abs > t & all.no.abs$prop.aye.pp.no.abs > t & all.no.abs$prop.aye.pr.no.abs > t, -ncol(all.no.abs)], 1)
# Left together (PT), against right-center (PMDB, PSDB, DEM, PR, PP)
round (all.no.abs[all.no.abs$prop.aye.pt.no.abs > t & all.no.abs$prop.aye.pmdb.no.abs < (1-t) & all.no.abs$prop.aye.pp.no.abs < (1-t) & all.no.abs$prop.aye.dem.no.abs < (1-t) & all.no.abs$prop.aye.psdb.no.abs < (1-t) & all.no.abs$prop.aye.pr.no.abs < (1-t), -ncol(all.no.abs)], 1)

round (all.no.abs[all.no.abs$prop.aye.pt.no.abs < (1-t) & all.no.abs$prop.aye.pmdb.no.abs > t & all.no.abs$prop.aye.pp.no.abs > t & all.no.abs$prop.aye.dem.no.abs > t & all.no.abs$prop.aye.psdb.no.abs > t & all.no.abs$prop.aye.pr.no.abs > t, -ncol(all.no.abs)], 1)

# Government in favor, opposition against
valid.vote[all.no.abs$V1[all.no.abs$prop.aye.pmdb.no.abs > t & all.no.abs$prop.aye.pt.no.abs > t & all.no.abs$prop.aye.pp.no.abs > t & all.no.abs$prop.aye.dem.no.abs < (1-t) & all.no.abs$prop.aye.psdb.no.abs < (1-t) & all.no.abs$prop.aye.pr.no.abs > t]]

# Government against, opposition in favor
valid.vote[all.no.abs$V1[all.no.abs$prop.aye.pmdb.no.abs < (1-t) & all.no.abs$prop.aye.pt.no.abs < (1-t) & all.no.abs$prop.aye.pp.no.abs < (1-t) & all.no.abs$prop.aye.dem.no.abs > t & all.no.abs$prop.aye.psdb.no.abs > t & all.no.abs$prop.aye.pr.no.abs < (1-t)]]

# Left-center together (PT, PMDB, PSDB) in favor, right against (DEM, PR, PP)
valid.vote[all.no.abs$V1[all.no.abs$prop.aye.pmdb.no.abs > t & all.no.abs$prop.aye.pt.no.abs > t & all.no.abs$prop.aye.psdb.no.abs > t & all.no.abs$prop.aye.dem.no.abs < (1-t) & all.no.abs$prop.aye.pp.no.abs < (1-t) & all.no.abs$prop.aye.pr.no.abs > t]]
# Left-center together (PT, PMDB, PSDB) against, right in favor (DEM, PR, PP)
valid.vote[all.no.abs$V1[all.no.abs$prop.aye.pmdb.no.abs < (1-t) & all.no.abs$prop.aye.pt.no.abs < (1-t) & all.no.abs$prop.aye.psdb.no.abs < (1-t) & all.no.abs$prop.aye.dem.no.abs > t & all.no.abs$prop.aye.pp.no.abs > t & all.no.abs$prop.aye.pr.no.abs > t]]

# Left together (PT) in favor, right-center against (PMDB, PSDB, DEM, PR, PP)
valid.vote[all.no.abs$V1[all.no.abs$prop.aye.pt.no.abs > t & all.no.abs$prop.aye.pmdb.no.abs < (1-t) & all.no.abs$prop.aye.pp.no.abs < (1-t) & all.no.abs$prop.aye.dem.no.abs < (1-t) & all.no.abs$prop.aye.psdb.no.abs < (1-t) & all.no.abs$prop.aye.pr.no.abs < (1-t)]]

# Left together (PT, PSDB) in favor, right-center against (DEM, PR, PP)
valid.vote[all.no.abs$V1[all.no.abs$prop.aye.pt.no.abs > t &  all.no.abs$prop.aye.pp.no.abs < (1-t) & all.no.abs$prop.aye.dem.no.abs < (1-t) & all.no.abs$prop.aye.psdb.no.abs > (1-t) & all.no.abs$prop.aye.pr.no.abs < (1-t)]]


# Left together (PT) against, right-center in favor (PMDB, PSDB, DEM, PR, PP)
valid.vote[all.no.abs$V1[all.no.abs$prop.aye.pt.no.abs < (1-t) & all.no.abs$prop.aye.pmdb.no.abs > t & all.no.abs$prop.aye.pp.no.abs > t & all.no.abs$prop.aye.dem.no.abs > t & all.no.abs$prop.aye.psdb.no.abs > t & all.no.abs$prop.aye.pr.no.abs > t]]


# Left together (PT, PSDB) against, right-center in favor (DEM, PR, PP)
valid.vote[all.no.abs$V1[all.no.abs$prop.aye.pt.no.abs < (1-t) & all.no.abs$prop.aye.pp.no.abs > t & all.no.abs$prop.aye.dem.no.abs > t & all.no.abs$prop.aye.psdb.no.abs < t & all.no.abs$prop.aye.pr.no.abs > t]]


# Left together (PT) against, right-center in favor (PMDB, PSDB, DEM, PR, PP)
valid.vote[all.no.abs$V1[all.no.abs$prop.aye.pt.no.abs < (1-t) & all.no.abs$prop.aye.pmdb.no.abs > t & all.no.abs$prop.aye.psdb.no.abs < (1-t) &
                           all.no.abs$prop.aye.dem.no.abs < (1-t)]]

valid.vote[all.no.abs$V1[all.no.abs$prop.aye.pt.no.abs >t & all.no.abs$prop.aye.pmdb.no.abs < (1-t) & all.no.abs$prop.aye.psdb.no.abs >t &
                           all.no.abs$prop.aye.dem.no.abs >t]]



count.all <- count.dem <- count.pt <- count.psdb <- count.pmdb <- count.pr <- count.pp <- c()
for (i in 1:ncol(unityParty)) {
  count.all[i] <- length(unityParty[,i][!is.na(unityParty[,i])]) 
  count.dem[i] <- length(unityParty[,i][unityParty[,i]=="DEM" & !is.na(unityParty[,i])]) 
  count.pt[i] <- length(unityParty[,i][unityParty[,i]=="PT" & !is.na(unityParty[,i])]) 
  count.pr[i] <- length(unityParty[,i][unityParty[,i]=="PR" & !is.na(unityParty[,i])]) 
  count.pp[i] <- length(unityParty[,i][unityParty[,i]=="PP" & !is.na(unityParty[,i])]) 
  count.psdb[i] <- length(unityParty[,i][unityParty[,i]=="PSDB" & !is.na(unityParty[,i])]) 
  count.pmdb[i] <- length(unityParty[,i][unityParty[,i]=="PMDB" & !is.na(unityParty[,i])]) 
}

num.Aye.function <- function (tmp) {
  no.Aye <- length (tmp[tmp==1 & !is.na(tmp)])
  no.Nay <- length (tmp[tmp==0 & !is.na(tmp)])
  aye.tot <- no.Aye/(no.Aye+no.Nay)
  return (aye.tot)
}

PMDB.all.aye <- PT.all.aye <- PSDB.all.aye <- DEM.all.aye <- PR.all.aye <- PP.all.aye <- c()
for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="PMDB" & !is.na(unityParty[,i]),i]
  PMDB.all.aye[i] <- num.Aye.function(tmp)
}

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="PT" & !is.na(unityParty[,i]),i]
  PT.all.aye[i] <- num.Aye.function(tmp)
}

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="PSDB" & !is.na(unityParty[,i]),i]
  PSDB.all.aye[i] <- num.Aye.function(tmp)
}

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="DEM" & !is.na(unityParty[,i]),i]
  DEM.all.aye[i] <- num.Aye.function(tmp)
}

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="PR" & !is.na(unityParty[,i]),i]
  PR.all.aye[i] <- num.Aye.function(tmp)
}

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="PP" & !is.na(unityParty[,i]),i]
  PP.all.aye[i] <- num.Aye.function(tmp)
}

# Calculate average majority bloc
majorityBloc <- function (x, y, z) {
  tmpData <- cbind (x,y,z)
  validVotes <- rowSums (tmpData)
  largestShare <- apply (tmpData, 1, max)
  tmp <- largestShare/validVotes
  return (tmp)
}

dem.maj <- majorityBloc (x=aye.vote.dem, y=nay.vote.dem, z=obs.vote.dem)
pt.maj <- majorityBloc (x=aye.vote.pt, y=nay.vote.pt, z=obs.vote.pt)
pr.maj <- majorityBloc (x=aye.vote.pr, y=nay.vote.pr, z=obs.vote.pr)
pp.maj <- majorityBloc (x=aye.vote.pp, y=nay.vote.pp, z=obs.vote.pp)
psdb.maj <- majorityBloc (x=aye.vote.psdb, y=nay.vote.psdb, z=obs.vote.psdb)
pmdb.maj <- majorityBloc (x=aye.vote.pmdb, y=nay.vote.pmdb, z=obs.vote.pmdb)

round (mean (dem.maj, na.rm=T), 2)
round (mean (pt.maj, na.rm=T), 2)
round (mean (pr.maj, na.rm=T), 2)
round (mean (pp.maj, na.rm=T), 2)
round (mean (psdb.maj, na.rm=T), 2)
round (mean (pmdb.maj, na.rm=T), 2)

round (min (dem.maj, na.rm=T), 2)
round (min (pt.maj, na.rm=T), 2)
round (min (pr.maj, na.rm=T), 2)
round (min (pp.maj, na.rm=T), 2)
round (min (psdb.maj, na.rm=T), 2)
round (min (pmdb.maj, na.rm=T), 2)

# Calculate average absence/abstention
round (mean (prop.abs.dem, na.rm=T), 2)
round (mean (prop.abs.pt, na.rm=T), 2)
round (mean (prop.abs.pr, na.rm=T), 2)
round (mean (prop.abs.pp, na.rm=T), 2)
round (mean (prop.abs.psdb, na.rm=T), 2)
round (mean (prop.abs.pmdb, na.rm=T), 2)


# Find Aye votes on items we use as anchors
anchor.votes <- c("Voto.2011026","Voto.2011075","Voto.2013077","Voto.2013078")
round (PT.all.aye[is.element(colnames(unityRollCall), anchor.votes)]*100, 1)
round (PMDB.all.aye[is.element(colnames(unityRollCall), anchor.votes)]*100, 1)
round (PP.all.aye[is.element(colnames(unityRollCall), anchor.votes)]*100, 1)
round (PR.all.aye[is.element(colnames(unityRollCall), anchor.votes)]*100, 1)
round (PSDB.all.aye[is.element(colnames(unityRollCall), anchor.votes)]*100, 1)
round (DEM.all.aye[is.element(colnames(unityRollCall), anchor.votes)]*100, 1)

unity.index.function <- function (tmp) {
  no.Aye <- length (tmp[tmp==1 & !is.na(tmp)])
  no.Nay <- length (tmp[tmp==0 & !is.na(tmp)])
  no.Abs <- length (tmp[is.na(tmp)])
  rice.wo.abs <- abs (no.Aye-no.Nay) / (no.Aye+no.Nay)
  rice.w.abs  <- abs (no.Aye-no.Nay) / (no.Aye+no.Nay+no.Abs)
  alt.rice.wo.abs <- max (no.Aye, no.Nay) / (no.Aye+no.Nay)
  alt.rice.w.abs  <- max (no.Aye, no.Nay) / (no.Aye+no.Nay+no.Abs)
  agreement.index <- (max (no.Aye, no.Nay, no.Abs) - 0.5*(length(tmp) - max (no.Aye, no.Nay, no.Abs))) / length(tmp)   
  all.unity.indexes <- c (rice.wo.abs, rice.w.abs, alt.rice.wo.abs, alt.rice.w.abs, agreement.index)
  return (all.unity.indexes)
}

# Build party-specific, vote-specific unity indices
PMDB.unity.index <- matrix (NA, nrow=5, ncol=ncol(unityRollCall))
PT.unity.index <- PSDB.unity.index <- DEM.unity.index <- PR.unity.index <- PP.unity.index <- PMDB.unity.index

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="PMDB" & !is.na(unityParty[,i]),i]
  PMDB.unity.index[,i] <- unity.index.function(tmp)
}

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="PT" & !is.na(unityParty[,i]),i]
  PT.unity.index[,i] <- unity.index.function(tmp)
}

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="PSDB" & !is.na(unityParty[,i]),i]
  PSDB.unity.index[,i] <- unity.index.function(tmp)
}

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="DEM" & !is.na(unityParty[,i]),i]
  DEM.unity.index[,i] <- unity.index.function(tmp)
}

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="PR" & !is.na(unityParty[,i]),i]
  PR.unity.index[,i] <- unity.index.function(tmp)
}

for (i in 1:ncol(unityRollCall)){
  tmp <- unityRollCall[unityParty[,i]=="PP" & !is.na(unityParty[,i]),i]
  PP.unity.index[,i] <- unity.index.function(tmp)
}

rownames (PMDB.unity.index) <- c("rice.wo.abs","rice.w.abs","alt.rice.wo.abs","alt.rice.w.abs","agreement.index")
rownames (PT.unity.index) <- rownames (PSDB.unity.index) <- rownames (DEM.unity.index) <- rownames (PP.unity.index) <- rownames (PR.unity.index)

#### Load party pressure full chain
#### The following line loads results from the Positive IRT pressure model
load ("brazil_ptirt_mcmc.Rda")

X.e   <- rbind (chainsPTIRT[[1]][,grep("X",  colnames (chainsPTIRT[[1]]))],
                chainsPTIRT[[2]][,grep("X",  colnames (chainsPTIRT[[2]]))])

rm (chainsPTIRT)

X.pt    <- X.e[,grep(",1]", colnames(X.e))]
X.pmdb  <- X.e[,grep(",2]", colnames(X.e))]
X.psdb  <- X.e[,grep(",3]", colnames(X.e))]
X.dem   <- X.e[,grep(",4]", colnames(X.e))]
X.pp    <- X.e[,grep(",5]", colnames(X.e))]
X.pr    <- X.e[,grep(",6]", colnames(X.e))]


# Uncomment next line to save file
# save.image("prepared_roll_call_data_brazil.Rdata")


